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Traffic jam in a fluid model of traffic flows proposed by Kerner and Konhauser (B. S. Kerner and P. Konhauser, 
Pliys. Rev. E 52 (1995), 5574.) is analyzed. An analytic scaling solution is presented near the critical point of the 
hetero-clinic bifurcation. The validity of the solution has been confirmed from the comparison with the simulation of 
the model. 
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Cooperative behavior in dissipativc systems composed 
of discrete elements has attracted much interest of physi- 
cists. Traffic flows have been studied extensively (f; 2; 
3; 4) as one of such systems. They can be classified 
into discrete models and continuous models. The for- 
mer includes the cellular automaton model (5) and the 
optimal velocity model(6), and the latter dose the fluid 
models. Recently Kerner and Konhauser (7) have pro- 
posed a fluid model of traffic flows which is regarded as a 
standard one. Since these models are aimed to describe 
common phenomena such as traffic-jam formation, they 
should have universal properties as stressed by Hayakawa 
and Nakanishi (8). The aim of this letter is to confirm 
the universality of traffic-jam formation. For this pur- 
pose we quantitatively compare the results for the fluid 
model with those for a discrete model (9). 

Komatsu and Sasa (10) have revealed the mechanism 
of jam formation in the optimal velocity model proposed 
by Bando et al. (6) with the aid of the perturbational 
treatment of solitons. Hayakawa and Nakanishi (9; 8) 
have developed the analysis of Komatsu and Sasa (10), 
and analyzed in detail a generalized optimal velocity 
model 

Xn = a\W{Xn+\ - Xn) V{Xn - Xn-\) - , (1) 

where a;„ and a are the positions of n-th car and the sen- 
sitivity, respectively. The driver care of not 
only the forward distance Xn+i — Xn but also the back- 
ward distance x„ — Xn-i- The optimal velocity function 
W is a monotonic increasing function, and V is a mono- 
tonic decreasing function. Hayakawa and Nakanishi (9; 
8) have obtained a scaled asymmetric kink solution of 
this model. On the other hand, Komatsu (11) has shown 
from his simulation that the fluid model by Kerner and 
Konhauser (7) has also a scaled asymmetric kink solu- 
tion as that in eq. (1). Here, we will clarify the rela- 
tion between the model by Kerner and Konhauser and 
eq. (1), based on the method developed by Hayakawa 
and Nakanishi (9; ?). As a result, we will demonstrate 
that the behavior of the fluid model dose not contain es- 
sential differences from that of the discrete models. In 
addition, the validity of our analysis will be confirmed 
from the comparison with our simulation of the model. 

Kerner and Konhauser proposed a one-dimensional 
compressible fiuid model for traffic flow in a highway 

( dt(j) = -d;,(f)v, 

If, , U{<p)-v T a . (2) 

1 dtv = -vd^v H — -r9z<P + -rd^v, ^ 

I. T (p (p 

where (/> is the field of car density, and v is the field of car 
velocity. The term {U{4)) — w)/t with the relaxation time 
T represents the relaxation process, where the optimal ve- 
locity field U corresponds to WV in eq. (1). The pressure 
term with the strength T and the viscous term with the 
effective viscosity /i stabilize the solution of the model, 
which is the effects implicitly included in the model (1). 
The form of the optimal velocity U{(j)) can be estimated 



from the relation between density and velocity of cars 
observed in real traffic ffows in highways (1). Its explicit 
form is assumed to be 

U{(p) = 2.52305 [tanh((l - 0.25)/0.12) 

-tanh(((;i-0.25)/0.12)], (3) 

as the paper by Kerner and Konhauser (7). 

First, we consider the linear stability of a uniform ffow, 
{4>,U = U{4>)), in the ffuid model (2). The linearized 
equation of the model (2) in the neighborhood of this 
uniform state is 

{dt4> = -(t>dzv, 

where (j) = 4'~'t>^''^ = v — f7oei?!oc(Band U' = dU/d(j)\^^^. 
The Fourier transform of eq. (4) satisfies 

where fc is a wave number. The growth rate a in the 
solution e*^* of this equation is the solution of 

^^ + (}L + iikAa+Te+i^k = o. (6) 

\T (j) J T 

This equation has two solutions. The condition for in- 
stability of the uniform state is Re{a+) > where a+ is 
one of the solution of eq. (6), 

Prom eq. (7) we find that A; = is the neutral mode and 
the mode fc — > is the most unstable. Thus, when 

^^U'^ > T (8) 

is satisfied, the uniform flow becomes unstable owing to 
emergence of the negative diffusion constant. 

For later convenience, we explicitly write the expansion 
of (7+ around k = 0. 

a+ = -i^U'k + t{-T + 4>'^U'^)k^ 

+ irU'in - 2tT4> + 2T4?U'^)k^ 

+ 6tT^^U'^ - 5r^^i7'*) k" + . . . . (9) 

Here, we notice that r and /x are not included in the 
instability condition (8). Then we assume r = 1 and 
^ = 1 by a proper transformation of the scales of the 
space, time and field quantities. 
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Next, wc carry out the weakly nonlinear analysis. 
Firstly, let us reduce the model (2) in the vicinity of 
neutral curve (8) with the aid of the reductive perturba- 
tion method (12). From eq. (9), the terms in proportion 
to k'^ and k'^, must be balance. Furthermore, since the 
term ik is eliminable by the Galilei transformation, the 
term ik^ is the dominant in the dispersion relation, which 
should be balance with the lowest order nonlinear term. 
If the lowest nonlinear term is quadratic, the model is 
reduced to the Korteweg-de Vries (K-dV) equation (13; 
14; 15; 16; 17). The K-dV equation, however, has only 
pulse solutions, so it is imsuitable for a description of the 
traffic jam formation. Therefore, we should choose the 
critical point which is the cross point of the neutral curve 
T = (l>'^U'{(t)f and (t>U"{(t>)+2U'{(t)) = where the cubic 
nonlinear term becomes dominant. The explicit critical 
point when we adopt eq. (3) is given by 



(0c, Tc) = (0.300704,28.2553). 



(10) 



Now let us look for a steady propagating solution in 
the neighborhood of this critical point. We introduce the 
scaled variable x 



(16) 



we obtain the following pair of asymmetric kink and anti- 
kink solutions in the lowest order. 



4'^\x)=ta.nhi0±xy, g± = ^^^i^'^ -. (17) 



Let us assume that the system satisfies the periodic 
boundary condition. Thus the solution should be modi- 
fied as 



Mx) ^ ^l^o^^\x - a;+) - 1 + iJo^->{x - X-). (18) 

This solution shows two interfaces at the positions x = 

x±, and approximately satisfies the periodic boundary 
condition, though we notice that eq. (18) is only an ap- 
proximate solution in the lowest order equation of (14). 
Substituting (18) into (14), we obtain 



(19) 



where 



(eiz - Vrt) - ceh^ 



(11) 



where e = y' {Tc — T)/Tc, A = —U', and c is the positive 
free parameter which will be determined from the pertur- 
bation analysis. Since the lowest nonlinear term becomes 
cubic, we assume following expansions: 



• = (l)c + e \/ -7p + (1)3 + e^04 + • • • , 
U{(pc) + evi + €^ V2 + V3 + Vi . . . , 



(12) 
(13) 



where B = <j)cU'^^^ /6 - U'/<j)c and -0 includes 0(e). We 
substitute these expression into eq. (2) and collect the 
terms in each power of e. From the expression to the 
fifth order of e, we obtain 



where 



^{i;^ - 1) + Pd,^^] = ed,M[^], (14) 



M[tp] = sTc 



c 



+ 



2ct 



-xipo — cct^ 



(15) 



where we assume c is 0(e), /3 = C/VAB, p23 = 

D/^A^, P32 = E/VAB^, Pii = F/Vb^, P14 = 

G/VA^. C = (f>eU'^, D = -2(j)lU'^ - U'I4>^, E = 

-<PlU'U<-^y3 + C/'2, F = -</>ef/W/24 - f/^^Ve, G = 
2(j)eU'^, and = </)2c/'2. 
Here, assuming the expansion 



+ 2(3dli,o + ^fidx'4>odx + 2/3Vo5^. 



(20) 



To obtain a regular behavior of perturbation in 0(e) 
the perturbed solution should be orthogonal to the zero 
eigenfunction (9). This solvability condition is given by 

{^o,dxM[ii;Q\)= \im [ dx^odxM[^o]=0, (21) 

where L is the system size, and is a zero eigenfunction 
which satisfies 

£t^-o = 0; = -dl-dx + 3i^ldx + 2pi;odl. (22) 
Equation (22) can be rewritten as 

C^<^o{x) =0; a = -dl-l + 3^Pl + 2Pi;odx, (23) 

where $o = c^x^o- A pair of special solutions of eq. (23) 
is given by 



^'^\x) = {sec]i{e±x)) 



(24) 



and the corresponding solutions of (22) can be expressed 
by 



^o^^Hx) 



dx'{sech{e±x')) 



(25) 



Since \l/o should satisfy the periodic boundary condition 
as well as ipo, we adopt the superposition of "^q^^ as 

*o(a;) = 'Ho'-^^x -x+)-l + ^o''~\x - X-). (26) 
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The constants a± in cq. (25) are determined to satisfy 
\l/o(±oo) = — 1. Thus, we choose 



a± = 



20± 



(27) 



where 



/oo 
dx{sechx)^ 
-OO 

. r(^ 



/ej.+2n 



(28) 



with the Gamma function r(a;). 

Now, we rewrite the solvabihty condition (21) by the 
integration by parts 



(29) 



From (15) it is obvious that terms in M[V'o] have no con- 
tribution in the left hand side of (29) except for terms in 
proportion to p^i and x-ijj. Notice that the contribution 
from the term in proportion to tip vanishes because of 
its symmetry. In addition, since we adopt the periodic 
boundary condition, as in eqs. (18) and (26), the contri- 
bution from the term p4i is canceled. Therefore the left 
hand side of (29) is reduced to 



[*oM[V;o]]; 



-L. 



C2 



(30) 



After computation of the right hand side of eq. (29), we 
obtain the time evolution equation of c, 



[L- {e+-e-)]c = 4pc^ 



1 



c 

c+ 



1-^ 

c_ 



where we use 

Here c± are given by 



2ne% + 1 
(2n+ 1)^1 + 1' 



(31) 



(32) 




(33) 



which are expected to be the stable values of c where 
a kink or an anti-kink exists under the free bound- 
ary condition. Equation (31) represents that the re- 
laxation time for a steady state is divergent in propor- 
tion to the system size L in the vicinity of the conver- 
gent value c* of e. From eqs. (3) and (10), each con- 
stant is determined as /3 = 2.01476, 9+ = 2.23815, 9- = 
-0.223398, P23 = 5.39424, pz2 = 1.92455, pn = 
0.299797, pi4 = 2.52857, pi2 = 6.72039. Therefore we 
obtain 



c+ = 2.67997, c_ 
and c* is given by 



2.62752. 



c = 



_ c+c-[e+{el + 1) - e-{9l + 1)] 



c-9+{9l + l)-c+9-{9% + l) 
Thus the explicit value of c* is 
c* = 2.66066. 



(34) 



(35) 



(36) 



To check the validity of our analysis we perform the 
numerical simulation of the model (2) in the vicinity of 
the critical point under the periodic boundary condition. 
We adopt the classical fourth-order Runge-Kutta scheme 
to time and the Eular scheme to space in our simulation. 
The analytic solution in eqs. (17), (18), and (36) was 
adopted as the initial condition for the fast convergence. 
We also have checked that the tendency to converged to 
c* from the random initial condition. Taking into account 
the scaling properties we perform the simulation for each 
parameter, e = 1/2,1/4,1/8, and 1/16, imtil (p relaxes 
to a steady propagating state. Here we change the space 
mesh size according to the divergence of the characteristic 
length in proportion to e ^ 0. 

Our results are plotted in Figs. 1 and 2. In Fig. 1, 
points represent the maximum and minimum values of 
car density in each parameter. The solid curve represents 
the theoretical coexistence curve which is obtained from 
the above perturbation analysis as 



T = Tr 1- 



c* 



(37) 



We emphasize that the deviation between simulation and 
theory is only 1.0% at e = 1/16. It is interesting that we 
can see one of the branches in the linearly unstable re- 
gion. Figure 2 displays the scaled density of cars as the 
function of the scaled position of cars. The solid line is 
theoretical one in eqs. (17), (18), and (36). The result of 
our simulation is asymptotically identical to our theoret- 
ical curve as e ^ 0, and reproduce an asymmetric pair 
of kink and anti-kink. As we have seen in this letter, our 
results of the simulation agree well with our theoretical 
analysis. 

In conclusion, our results in the fluid model in eq. (2) 
can be summarized as follows: 
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1. The scaling solution exists in the vicinity of the 
critical point. 

2. There is an asymmetry of the solution. 

3. Quantitative results can be reproduced by the per- 
turbation analysis. 

These results correspond to those in the analysis of dis- 
crete model (1)(9; 8). Thus we suggests that the results 
obtained here are universal ones which are independent 
of choice of a specific model. 
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FIG. 1. The theoretical coexistence curve (sohd hne) in 
eq. (31) and the neutral curve (dotted line) in eq. (8) in (0, T) 

plane. The data points (cross) arc obtained from the maxi- 
mum and minimum values of car density </) at given T. 

FIG. 2. The theoretical curve in eqs. (17), (18), and (36) 
and data of the scaled density obtained by our simulation. 
The horizontal axis is displays the car position. The system 
size is normalized to 1, and positions of the interfaces are 
set to 0.25 and 0.75 which are the fitting parameters of the 
theoretical expression. Each data is obtained from the result 
for e = 1/2, 1/4, 1/8, and 1/16. 
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